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ABSTRACT 


Modeling  and  simulation  of  crack  initiation  and 
propagation  in  solid  rocket  propellant  materials  were 
conducted  using  both  the  macromechanics  approach  and  the 
micro /macromechanics  approach.  Due  to  their  composition, 
the  solid  rocket  propellant  can  be  construed  as  particle 
reinforced  composites. 

The  macromechanics  approach  entailed  a  numerical 
simulation  of  a  finite  element  model  to  predict  the  crack 
behavior  based  on  the  damage  initiation,  growth,  and  local 
saturation.  Its  results  were  then  compared  to  the 
experimental  data.  In  the  simulation,  it  was  assumed  that  a 
crack  forms  when  a  damage  is  saturated  in  a  localized  zone. 
The  results  from  the  simulation  were  quite  comparable  to  the 
experimental  results,  validating  the  method  of  predicting 
crack  initiation,  growth,  and  arrest  using  the  concept  of 
damage  growth  and  saturation. 

The  second  approach  involved  using  a  simplified 
micromechanical  model  and  the  damage  mechanics  being  applied 
at  the  micromechanics  level  and  the  finite  element  analysis 
being  done  subsequently  at  the  macromechanics  level.  In 
using  this  approach,  the  damage  modes  such  as  matrix 
cracking,  interface  debonding  and  particle  cracking  were 
explainable  in  an  explicit,  fundamental  manner.  Several 
simulations  were  conducted  using  this  approach  including  the 
cases  of  non-uniform  particle  distribution.  The  predicted 
results  comapred  well  with  the  experimental  data. 
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I. 


INTRODUCTION 


Because  of  their  light  weight  and  stability,  in  the 
last  few  decades,  solid  propellants  have  replaced  liquid 
fuel  as  the  primary  means  of  projecting  rockets  and 
missiles.  Most  missiles  in  use  today  such  as  the  Tomahawk 
and  the  Harpoon  employ  solid  propellants.  Missiles  like 
these,  because  of  the  nature  of  their  use,  require  high 
level  of  accuracy  and  reliability.  This  means  all  systems 
onboard  the  missile  are  required  to  operate  exactly  to  he 
designed  specifications  including  the  propulsion  system. 
Solid  propellants,  despite  their  advantages,  have  few 
problem  areas.  One  of  these  areas  is  the  initiation  of 
damage  and  cracks  prior  to  and  during  operation.  Microflaws 
over  time  tend  to  soften  the  mechanical  behavior  and 
increase  the  potential  for  the  formation  of  a  macroscopic 
crack.  When  this  happens,  the  increased  burning  surface 
area  as  a  result  of  the  crack  maybe  enough  to  over¬ 
pressurize  and  thereby  cause  failure  in  the  motor  casing 
[Cornwell  &  Schapery  (1975)].  Even  if  failure  in  the 
casing  may  not  occur,  the  formation  of  the  crack  can  cause 
uneven  burning  of  the  propellant  leading  to  decrease  in 
range  of  the  missile. 

The  composition  of  solid  propellants  is  such  that  it 
can  be  construed  as  a  particle  reinforced  composite.  More 
specifically,  these  composites  are  made  up  of  a  mixture  of 
relatively  soft,  ductile  matrix  material,  reinforced  by 
randomly  scattered,  stiff  particles.  The  advantages  of 
these  composites  are  increase  in  elastic  modulus,  increased 
strength,  and/or  increased  toughness,  all  of  which  in 
varying  degrees.  [ Ravichandran  &  Liu  (1994)]  The  drawback 
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to  these  composites  is  that  unlike  the  increase  in  stiffness 
or  toughness,  the  benefit  of  the  increase  in  strength  is 
quite  limited.  [Johnson  (1987)]  As  a  result,  it  is 
critical  to  understand  and  to  predict  the  failure  mechanisms 
and  the  failure  modes  in  the  particle  reinforced  composites. 

To  understand  the  failure  modes  however,  the 
strengthening  factors  must  first  be  considered.  Extensive 
studies  by  Bretchet  (1990)  revealed  that  the  strength  of  the 
combined,  composite  material  depends  not  only  on  the 
inherent  material  properties  of  the  reinforcing  particles 
but  also  on  the  size  and  the  aspect  ratio  of  these 
particles.  They've  concluded  that  the  degree  of  influence 
of  these  factors  depend  on  the  Weibull  modulus  of  the 
particles  and  the  influence  of  these  factors  are  greater  in 
ductile  materials  as  opposed  to  brittle  materials.  Similar 
study  done  by  Shen  (1993)  in  the  effects  of  particle 
reinforcements  to  the  effective  elastic  modulus  of  the 
composite  material,  showed  that  the  cross-sectional  shape 
and  the  spatial  distribution  of  the  reinforcing  phase 
greatly  influenced  the  elastic . modulus.  More  specifically, 
it  was  found  that  in  testing  different  shapes  of  reinforcing 
agents  (of  aspect  ratio  approximately  equal  to  unity),  the 
effective  Young's  modulus  (at  a  fixed  concentration  of  the 
phase)  increased  in  the  following  order:  double  cone  — > 
sphere  — >  Octagon  — >  unit  cylinder.  [Shen,  Finot, 
Needleman,  &  Suresh  (1993)]  As  for  the  distribution,  it  was 
found  that  an  even  distribution  of  the  reinforcing  agent 
throughout  the  composite  material  gave  the  greatest  rise  to 
the  elastic  modulus. 
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The  objectives  for  this  research  were  to  develop  a 
finite  element  model  to  accurately  describe  the  experimental 
samples  of  particulate  composites  undergoing  uniaxial 
loading  and  to  develop  simulation  techniques  to  investigate 
damage,  crack  initiation,  and  growth  in  particulate 
composites.  By  providing  an  accurate  method  of  predicting 
the  crack  initiation  and  propagation,  the  need  for  costly 
and  time  consuming  experiments  could  be  decreased.  Also, 
parameters  such  as  the  loading  condition  and  the  material 
property  can  be  changed  easily  for  a  computer  simulation 
than  for  an  actual  experiment.  Furthermore,  these 
information  help  engineers  and  scientists  understand  the 
mechanisms  of  fracture  in  particle  composites. 
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II.  BACKGROUND  LITERATURE  REVIEW 


In  this  study,  two  methods  are  used  to  model  the  crack 
initiation  and  growth.  Although  the  specific  procedure  and 
the  methods  used  in  this  study  are  unique,  there  have  been 
similar  studies  done  in  the  past  that  agree  in  principle  and 
approach.  In  determining  the  effective  properties  of 
heterogeneous  materials,  Hashin  (1962)  utilizes  variation 
theorems  of  the  theory  of  elasticity  to  achieve  upper  and 
lower  bounds  for  the  Young's  modulus.  In  his  work,  Hashin 
summarizes  that  in  determining  the  elastic  moduli  of 
heterogeneous  materials,  they  can  be  described  by  an 
elastic,  isotropic  matrix  whose  moduli  are  known.  These 
isotropic  particles  are  embedded  with  known  homogenous 
elastic  properties.  Assuming  also,  that  the  volume  fraction 
of  the  particles  are  evenly  distributed  and  as  a  result,  be 
regarded  as  quasi-homogenous  and  quasi-isotropic,  the 
problem  becomes  finding  the  elastic  moduli  of  this  combined 
material.  Specifically,  his  method  derives  the  effective 
moduli  involving  only  the  stresses  or  strains  inside  the 
inclusions.  This  is  done  by  considering  the  change  in 
strain  energy  in  a  loaded  homogenous  body  due  to  the 
inclusion  of  nonhomoeneities .  He  begins  his  derivation  by 
applying  the  Hooke's  law  to  the  strain  energy  equation  for 
an  elastic,  homogenous,  isotropic  material  and  then 
splitting  the  stress  and  strain  tensors  into  isotropic  and 
deviatoric  parts ,  as  shown  in  Equation  ( 1 ) . 
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where 


a^°'> 

O  kk  f 


^ij  ij  2  ^  ^  'J 


and  K„  is  the  bulk  modulus,  G„  is  the  shear  modulus,  (j y  is 
the  stress  tensor,  and  by  is  the  Kronecker  delta.  Then, 
after  including  N  particles  into  the  material  and  applying 
surface  tractions  (tensile  load  in  the  vertical  direction  & 
compressive  load  in  the  horizontal  direction) ,  the  strain 
energy  stored  in  the  composite  body  can  be  expressed  by 

=Uo+bU^°^  =Uo+'^bU„^°^  (2) 

n=l 


where 


where 


o(0)  _p(0) 

^ij  ~^kk 

AO)  _p(0)  _ig(0)g 

^ij  -^ij  3  ^  °  ij 


The  notations  with  subscript  of  "0"  refers  to  the  same  given 
body  but  with  no  inclusions  and  under  the  same  load.  Zy 

denotes  the  strain  tensor,  Vn  denotes  the  volume  of  the  nth 
inclusion,  and  the  subscript  "p"  denotes  the  particle 
material  and  the  corresponding  stresses  and  strains  are 
written  without  superscripts .  By  the  above  equations 
therefore,  if  the  stresses  and  strains  in  the  inclusions  are 

known,  then  can  be  determined.  Now,  assuming  the 
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tractions  loaded  on  the  boundary  surface  is  derived  from  a 
uniform  stress  system,  can  then  be  expressed  by 

Equation  (2)  where  the  stresses  are  constant  and  the  moduli 
and  G„  can  be  replaced  by  the  effective  moduli  K'  and 
G* .  Therefore, 


(  (0)' 

a 


^  (0)^  (0)^ 
2G* 


(4) 


where  the  explanation  for  various  quantities  are  as 
previously  given.  Knowing  then,  the  above  equation  can 
then  be  expressed  in  terms  of  the  effective  moduli-  [Hashin 
(1962)] 

The  damage  modes  of  a  particle  reinforced  composite  are 
particle-matrix  debonding  called  dewetting,  matrix  damage 
such  as  cavity  formation,  and  particle  cracking.  If  the 
difference  in  strength  is  significantly  large  between  the 
matrix  material  and  the  reinforcing  particles,  the  mode  of 
damage  or  failure  is  usually  dewetting  or  cavity  formation 
in  the  matrix  or  both.  Studies  done  by  Schapery ( 1987 ) , 
Anderson  and  Farris  (1988),  and  Ravichandran  and  Liu  (1995) 
have  shown  that  these  damages  cause  volume  dilatation 
resulting  in  a  nonlinear  stress-strain  behavior.  In  other 
words ,  they  investigated  the  nonlinear  constitutive  response 
of  the  damaged  particulate  composite  with  the  change  in 
volume  dilatation  as  a  weighing  factor.  Micromechanical 
models  were  considered  to  describe  the  nonlinear  stress- 
strain  behavior  with  damage  evolution  in  particle  reinforced 
composites.  These  studies  also  included  the  failure  modes. 
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dewetting  and  matrix  cracking.  Graphic  representation  of 
these  failure  modes  are  included  in  Figure  1. 
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The  application  of  micromechanics  in  predicting  the 
effective  properties  and  failure  characteristics  of 
composite  materials  is  not  unique  to  the  present  study. 
There  in  fact,  exists  several  well-known  models  that  have 
been  developed,  and  Christensen  brings  these  models  together 
in  a  critical  evaluation  in  regard  to  their  range  of 
applicability  and  their  accuracy.  The  three  micromechanics 
models  which  he  considers  are  the  Differential  Method 
[Roscoe  (1952)],  the  Generalized  Self  Consistent  Method 
[Christensen  and  Lo  (1979)],  and  the  Mori-Tanaka  Method 
[Benveniste  (1987)].  His  evaluation  of  the  above  methods  is 
on  a  material  with  poly-dispersed  spherical  inclusions  whose 
limiting  case  is  specified  by  those  orientations  that  allow 
complete  packing  with  v.f.  — >  1.  There  are  several 
advantages  for  considering  models  that  permit  v.f.  — >  1. 
First  is  that  they  are  well  suited  for  theoretical 
processing.  As  a  result,  they  are  widely  used  in 
theoretical  studies.  Second  is  that  they  model  the  actual, 
poly-disperse  particle  composites  quite  accurately,  and  the 
third  is  that  they  are  an  excellent  approximations  to  the 
mono-disperse  case  for  concentrations  that  are  not  too 
large.  Christensen  begins  his  evaluation  with  the 
Differential  Method.  The  concept  behind  this  method  is  to 
view  the  whole  composite  as  a  sequence  of  dilute  suspensions 
and  applying  the  below  equations  for  the  effective  shear  and 
bulk  moduli  (of  an  isotropic  composite  with  non-interacting 
spherical  inclusions)  as  each  inclusion  is  added  one  at  a 
time. 
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where  subscripts  i  and  m  denote  inclusion  and  matrix,  v  is 
the  Poisson's  ratio,  and  k  denote  effective  shear  and  bulk 
modulus  respectively,  and  c  denotes  the  volume  fraction  of 
the  particles.  Therefore,  after  calculating  pi  and  k  with 
one  inclusion,  that  suspension  is  considered  as  a  homogenous 
material  with  those  properties,  to  which  the  second 
inclusion  is  added,  and  the  calculation  is  done  again.  The 
process  continues  until  the  condition  of  full  packing  of  the 
inclusion  is  reached.  The  final  result  of  this  differential 
process  is  the  set  of  governing  equations: 
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where  effective  Poisson's  ratio  v  is  equal  to 
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The  boundary  conditions  for  the  above  equations  are 

at  c=0:  n  =  (10) 

k  =  km 


at  c  =  l:  ^  =  ^,- 
k  =  kj 

The  second  method  Christensen  describes  is  the  Generalized 
Self  Consistent  Method.  The  concept  behind  this  method  is 
to  solve  the  equations  derived  from  the  model  where  a 
spherical  inclusion  is  placed  concentrically  inside  a  matrix 
material  which  in  turn  are  embedded  in  an  effective  medium 
with  the  properties  of  interest,  and  k.  The  model  is 
mathematically  represented  by  a  quadratic  equation, 


+  2B\ 


+  c=o 


(11) 


where  the  coefficients  A,  B,  and  C  are  equal  to  complex  and 
rather  tedious  equations  involving  Poisson's  ratio  of  the 
matrix  and  shear  moduli  of  the  inclusions  and  the  matrix. 
The  final  result  of  the  above  quadratic  equation  is  a  closed 
form  solution  for  the  effective  modulus. 
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The  final  micromechanics  analysis  Christensen  presents  is 
the  Mori-Tanaka  Method.  The  Mori-Tanaka  Method  is 
significantly  different  from  the  first  two  models  in  that 
where  as  the  first  two  methods  concede  physical 
descriptions,  the  Mori-Tanaka  Method  relies  on  mathematical 
manipulation  of  the  field  variables  along  with  special 
concepts  of  eigenstrain  and  backstress .  For  the  development 
of  this  method,  Christensen  cites  Benveniste  (1987)  who 
developed  a  simplified  derivation  of  the  method.  (However, 
even  the  simplified  method  is  fairly  extensive.)  Since  the 
derivation  is  quite  involved,  just  the  final  result  will  be 
presented.  The  manipulation  of  Mori-Tanaka  Method  yields. 


k  =  k^+cik,-k„) 
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l  +  (l-c) 
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(13) 


Having  presented  the  details  of  the  three  micromechanics 
models,  Christensen  then  evaluates  the  three  methods  in 
terms  of  their  compatibility  to  the  experimental  results. 
The  conclusion  he  makes  based  on  the  evaluation  is  that  the 
Generalized  Self  Consistent  Method  is  the  most  accurate  of 
the  three  methods  in  determining  the  effective  properties  of 
the  composite  material,  especially  in  the  limiting  case  of  a 
concentrated,  poly-disperse  suspension  of  rigid  particles. 
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The  other  two  models  gave  inaccurate  and  unrealistic 
behavior  in  the  limiting  case.  [Christensen  (1989)] 


III.  PRE-ANALYSIS  PROCEDURE  AND  MODELING 


A.  PROCEDURE 

Initially,  a  2-D  Finite  Element  Model  was  developed 
using  the  IDEAS  program  produced  by  SDRC.  Once  the  models 
were  generated,  the  mesh  and  nodal  connectivity  information 
was  converted  to  an  "ASCII"  universal  file  and  extracted 
from  IDEAS.  This  information  was  read  into  a  locally 
developed,  general  purpose  FORTRAN  finite  element  analysis 
program  which  could  be  set  up  in  various  ways  to  provide 
desired  information  such  as  effective  stress,  displacement, 
principle  stress,  micro-stresses  etc.  In  executing  this 
step,  however,  some  modifications  had  to  be  made  to  the 
input  data.  Specifically,  the  universal  output  file  from 
the  IDEAS  program  had  to  be  fitted  into  a  format  which  the 
FORTRAN  program  could  understand.  This  was  done  using  a 
simple  FORTRAN  code  which  discarded  unnecessary  information 
from  the  universal  file  and  place  the  remaining  data  in  a 
suitable  form.  Some  problems  that  were  encountered  during 
this  process  included;  (1)  the  local  nodal  mapping  which 
IDEAS  used  was  in  the  direction  opposite  to  the  FORTRAN  FEA 
program  and  (2)  the  rotation  constraints  had  to  be  left  at  a 
value  of  "0"  (i.e.  no  constraint)  for  the  FEA  program  to 
work.  Identifying  and  solving  these  compatibility  problems 
between  the  IDEAS  program  and  the  FORTRAN  FEA  program  posed 
significant  obstacles  in  the  progress  of  the  research.  Once 
the  mesh  information  was  successfully  entered  into  the 
FORTRAN  program,  the  rest  of  the  input  data  such  as  the 
boundary  conditions,  applied  loads,  and  material  properties 
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were  entered.  The  details  of  determining  these  input 
conditions  is  discussed  in  the  next  section. 

B.  FINITE  ELEMENT  MODELING 

The  driving  condition  for  the  FEM  models  in  terms  of 
their  physical  description  was  the  experimental  samples. 
The  dimensions  of  the  experimental  samples,  and  hence,  the 
FEM  models,  are  shown  in  Figures  2  and  3.  As  can  be  seen, 
two  samples  were  experimented.  One  was  a  3  in.  by  2.5  in. 
sample  with  a  .5  in.  diameter  hole,  and  the  other  was  a  3 
in.  by  2.5  in.  sample  with  a  .25  in.  diameter  hole.  These 
samples  were  placed  in  an  uniaxial  loading  machine  and  were 
loaded  to  .2  in/in  strain.  To  model  the  experiment,  a 
quarter  section  of  the  experimental  samples  was  used,  as  can 
be  seen  in  the  figures.  As  for  the  finite  element  mesh,  due 
to  the  constraints  in  the  computer  memory,  there  were 
limitations  in  number  of  nodes  and  elements  that  could  be 
used.  For  both  of  the  sample  models,  750  elements  and  806 
nodes  were  used.  As  was  found  later  during  the 

computational  phase  however,  the  degree  of  freedom  of  the 
mesh  was  still  too  large  for  the  computer  to  handle. 
Therefore,  further  manipulation  of  the  mesh  was  necessary. 
In  the  IDEAS  program,  there  exists  several  functions  which 
sequentialize  the  nodes  in  a  way  that  decreases  the 
bandwidth  of  the  global  matrix  used  in  the  calculations. 
After  several  different  variations  of  nodal  sequencing,  it 
was  found  that  sweeping  the  mesh  through  y,  then  x-axis  and 
subsequently  applying  the  Sloan  profile  reduction  scheme  was 
optimum.  By  doing  this,  the  bandwidth  decreased  from  58  to 
27  which  was  within  the  capability  of  the  computer.  The 
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Figure  2 .  Finite  Element  Model  of  Larger  Hole  Sample 
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next  step  in  the  modeling  process  was  to  apply  proper 
constrains  and  boundary  conditions.  Since  the  FEM  model 
represented  a  quarter  section  of  the  actual  experimental 
sample,  in  order  to  accurately  describe  the  experimental 
conditions,  the  left  vertical  edge  was  constrained  in  the  x- 
direction,  and  the  bottom  horizontal  edge  was  constrained  in 
the  y-direction,  with  the  left  corner  point  being 
constrained  in  both  x  and  y-directions .  The  load  was  placed 
uniformly  on  the  top  vertical  edge  of  the  model.  This  is 
shown  graphically  in  Figure  4 .  As  for  determining  the 
behavior  of  the  material,  two  conditions  were  considered, 
the  geometric  linear  and  the  nonlinear.  Initial  studies 
indicated  that  the  results  from  the  nonlinear  model  was  more 
realistic  in  terms  of  crack  size  for  a  given  load  and  the 
shape  of  the  damage  zone  near  a  stress  concentration  area. 
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Figure  4  -  FE  Model  with  Constraints  and  Loads  Shown 
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IV.  MACROMECHANICS  MODELING  AND  ANALYSIS 


A.  MODELING  AND  APPROACH 

In  order  to  simulate  the  experiment  conducted  for  the 
study  of  damage/crack  initiation  and  propagation  of  a 
particle  reinforced  composite  material,  the  finite  element 
method  was  used.  The  FEM  used  in  the  study  was  dictated  by 
the  nonlinear  elastic  behavior  of  the  composite  material. 
Figure  5  shows  a  typical  stress-strain  plot  of  a  particle 
reinforced  composite  material  having  undergone  a  uniaxial, 
tensile  test.  The  figure  also  shows  the  volume  dilatation 
of  the  material  as  a  function  of  strain.  The  figure  shows 
that  the  stress-strain  curve  begins  linearly,  and  as  the 
volume  dilatation  begins  to  increase,  it  starts  to  behave 
nonlinear ly.  After  the  continuous  nonlinear  behavior,  the 
stress-strain  curve  reaches  the  peak  point  and  the  failure 
of  the  specimen  occurs.  The  graph  also  shows  the  monotonic 
increase  in  volume  dilatation  until  the  peak  point  of  stress 
value.  This  increase  in  volume  dilatation  results  from  the 
formation  of  cavities  in  the  composite  material  during  its 
continuous  deformation.  Therefore,  the  nonlinear  behavior 
of  the  stress-strain  curve  can  be  attributed  to  these  micro¬ 
damages  in  the  material  such  as  formation  of  cavities .  In 
this  regard  therefore,  the  starting  point  on  the  stress- 
strain  curve  where  nonlinearity  begins  to  take  effect  is  the 
point  at  which  the  damage  or  crack  initiates.  On  the  same 
token,  the  peak  stress  point  on  the  curve  represents  the 
damage  saturation  state.  These  two  points  are  marked  in 
Figure  5 . 
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Stress-Strain  and  Dilatation  Curves 


DILATATION 


The  nonlinear  stress-strain  curve  is  used  in  the  finite 
element  analysis  in  order  to  model  the  damage  initiation, 
growth,  and  saturation  in  the  composite  material.  The  curve 
was  discretized  into  a  piecewise  linear  curve  which  was  then 
applied  to  the  finite  element  analysis .  The  actual  stress- 
strain  curve  used  in  the  numerical  analysis  is  shown  in 
Figure  6 .  In  extending  the  uniaxial  stress  information 
given  by  the  stress-strain  curve  to  the  multiaxial  state  of 
stress,  which  is  the  case  in  the  experiment,  the  von  Mises 
equivalent  stress  is  utilized  as  shown  below. 


a  +(^^2 +(^^3 


(14) 


When  the  equivalent  stress  reaches  the  peak  stress  value  at 
a  local  area,  the  area  is  considered  to  have  reached  the 
damage  saturation  point.  As  a  result,  when  this  occurs,  the 
damaged  material  cannot  sustain  any  additional  load. 
Therefore,  majority  of  the  load  is  transferred  away  from  the 
damaged  zone  to  the  adjacent,  undamaged  or  partially  damaged 
zones . 


B.  ANALYSIS  OF  THE  EXPERIMENTAL  DATA 

Figure  7  shows  the  X-ray  images  of  the  crack 
propagation  in  a  particulate  composite  specimen  under  mode  I 
condition,  uniaxial,  tensile  load  along  the  y-axis.  In  the 
image,  the  crack  is  represented  by  the  black  color  where  as 
the  specimen  irt  the  undamaged  areas  is  shown  in  color  white. 
The  dark  cluster  of  circular  spots  indicate  the  damage  zone. 
The  lighter  cluster  inside  the  cluster  which  is  most 
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(psi) 


clear  in  the  second  image  represent  a  more  damaged  zone.  As 
is  visible  in  Figure  7(a),  the  damage  zone  develops  and 
grows  for  the  most  part,  in  the  direction  of  the  crack.  As 
the  load  increases,  the  damage  zone  grows  and  spreads  in  the 
vertical  direction  as  seen  in  Figure  7(b).  Notice  the 
concentrated  damage  zone  just  ahead  of  the  crack  tip 
indicated  by  a  lighter  shade  as  the  crack  tip  becomes  blunt. 
Once  the  damage  becomes  saturated  in  front  of  the  blunt 
crack  tip,  the  crack  propagates  further  protruding  its  tip 
as  sharply  as  before  prior  to  the  blunting.  Therefore,  as 
it  is  seen  clearly  from  the  images,  when  the  damage  becomes 
saturated  in  a  local  region,  the  crack  forms  in  the  area. 
This  again  can  be  attributed  to  the  coalescing  of  the  micro¬ 
damages  in  the  damage  zone  and  forming  a  macro  size  crack. 
Hence,  in  the  analysis  of  the  FEM  model,  this  reasoning  was 
followed  in  determining  the  formation  of  the  crack.  In  the 
analysis,  the  particle  reinforced  composite  was  assumed  to 
behave  elastically.  In  other  words,  the  viscous  effect  was 
neglected,  and  the  strain  rate  effect  was  not  taken  into 
account.  Including  these  effects  would  be  the  logical  next 
step  in  the  follow  on  studies. 
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V. 


MICRO/MACROMECHANICS  MODELING  AND  ANALYSIS 


A.  SUMMARY  OF  METHOD 

The  developed  micro/macromechanical  method  utilizes  a 
micro-level  analysis  in  conjunction  with  a  macro-level 
analysis,  hence  the  name.  The  micro-level  analysis  includes 
computation  of  stresses  and  strains  in  the  constituent 
materials,  namely  the  reinforcing  particulates  and  the 
binding  matrix.  The  macro-level  analysis  involves 
structural  analysis  of  the  composite  as  a  whole.  As  will  be 
revealed,  the  two  analyses  were  performed  at  different 
scales  within  the  model  but  maintained  interaction  with  one 
another.  A  graphical  representation  of  the  interaction 
between  the  two  analyses  is  shown  in  Figure  8 .  The 
macromechanical  analysis  utilizes  the  finite  element  method 
for  structural  analysis  of  a  composite  in  a  way  in  which  the 
composite  structure  can  be  analyzed  generically.  In  other 
words,  the  analysis  was  conducted  so  that  the  solution  can 
be  in  multiple  forms  depending  on  the  information  needed. 
In  the  micromechanical  analysis,  a  simplified  unit-cell 
model  [Aboudi  (1987,1989);  Kwon  (1993):  Kwon  and  Berner 
(1994,1995)]  is  used.  Also,  on  the  unit-cell  model,  a 
damage  mechanics  [Kachanov  (9172,1980);  Krajcinovic  (1983); 
Talreja  (1985)]  is  applied.  The  simplified  unit-cell  model 
calculates  stresses  and  strains  at  the  constituent 
materials,  and  the  damage  mechanics  is  applied  to  the 
micromechanical  model  based  on  the  found,  micro-level 
stresses  and  strains.  One  advantage  of  this  approach  is 
that  the  damage  modes  and  mechanisms  can  be  described 
discriminately ,  in  terms  of  the  constituent  materials.  As  a 
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result,  the  damage  modes  can  be  specified  as  matrix 
cracking,  particle  cracking,  or  particle /matrix  debonding. 
For  the  sake  of  emphasizing  its  exactness  and  the  degree  to 
which  the  analysis  goes  in  depth,  it  is  of  interest  to 
describe  the  interaction  between  the  micromechanical  and 
macromechanical  analyses  in  detail.  The  micromechanical 
analysis  computes  the  smeared  composite  material  properties 
based  on  the  constituent  material  properties  and  their 
damage  states.  This  computation  is  conducted  at  every 
integration  points  of  the  macro-level,  finite  element  model. 
Therefore,  the  finite  element  analysis  can  be  undertaken  as 
usual  using  the  smeared  composite  material  properties  at  the 
given  integration  points .  After  the  macromechanical 
analysis  is  performed,  the  smeared  stresses  and  strains  at 
the  composite  level  are  then  applied  to  the  micromechanical 
model  so  that  micro-level  stresses  and  strains  can  be 
computed.  This  computation  is  performed  at  the  same 
integration  points  as  before.  Therefore,  as  a  result  of 
fixing  the  points  at  which  the  numerical  calculations  are 
conducted,  the  interaction  between  the  micro-level  analysis 
and  the  macro- level  analysis  occurs  at  the  integration 
points  of  the  finite  element  model.  And,  the 
macromechanical  model  represents  the  material  around  each 
integration  point.  These  domains  around  the  integration 
points  are  dependent  of  the  size  of  the  finite  element. 
Therefore,  the  damage  state  computed  from  the 
micromechanical  analysis  also  indicates  average  damage  of 
the  area  surrounding  each  integration  point. 

It  should  be  apparent  by  now  then,  that  by  using  the 
micromchanical  approach,  a  particle  reinforced  composite 
with  non-uniform  distribution  can  easily  be  modeled. 
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B.  UNIT-CELL  MICROMECHANICAL  MODEL 

The  microitiechanical  model  is  based  on  a  unit  cell 
structure.  A  unit  cell  consisting  of  a  particle  and  the 
surrounding  matrix  material  is  modeled  as  shown  in  Figure  9. 
The  figure  shows  a  two-dimensional  model.  This  is  an 
adoption  of  the  three-dimensional  model  which  was  studied  by 
Kwon  and  Berner.  (1994,1995)  Due  to  symmetry,  a  quarter 
section  of  the  full  unit— cell  model  is  considered,  and 
that's  what  is  shown  in  the  figure.  The  quarter  section  of 
the  unit-cell  consists  of  four  subcells  of  which  the  bottom 
left  cell  represents  the  particle  inclusion  and  the  rest  of 
the  surrounding  three  subcells  represent  the  matrix 
material.  By  this  model,  it  is  assumed  that  the  average 
deformation  of  the  unit  cell  represents  that  of  the 
particulate  composite  at  each  location. 

The  micromechanical  model  has  two  objectives:  One  is 
to  compute  the  smeared  composite  material  properties  from 
particle  and  matrix  material  properties  under  an  intact  or  a 
damaged  state.  The  other  is  to  describe  damage  modes  and 
states  at  the  specific,  constituent  material  level.  The 
development  of  the  micromechanical  model  is  described  below. 
As  stated  before,  this  two-dimensional  model  is  derived  from 
the  three-dimensional  model  which  was  developed  previously. 
Each  subcell  has  a  constitutive  equation  expressed  as 


where  i,j,k,l  =  1,2  and  a  =  a,b,c,d  Here,  1  and  2  denote 
coordinate  axes  while  a,b,c  and  d  indicate  subcells. 
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Subcell  a  represents  the  particle  inclusion  and  the  others 
represent  the  surrounding  matrix  material.  and  8^“  are 

the  average  stresses  and  strains  of  the  subcell  a.  The 
unit  cell  stresses  and  strains  are  obtained  from  the  volume 
average  of  the  subcell  stresses  and  strains.  That  is, 

5',=ro/  +  KV/  +  rcT,/+r‘'<T,/  (16) 

(17) 

where  O/j  and  Sy  are  the  unit  cell  stresses  and  strains 
which  eventually  represent  the  smeared  stresses  and  strains 
of  the  particulate  composite.  If  there  is  no  void  present, 
the  subcell  volume  fractions  can  be  computed  from  the 
particle  volume  fraction  bases  on  Figure  8. 

V"  =V^  (18) 

V'’  =V‘ 

v'' 


in  which  V’’  is  the  particle  volume  fraction. 

One  of  the  objectives  of  the  micromechanical  model  is 
to  find  the  relation  between  the  unit  cell  stresses  and 
strains  (i.e.  the  constitutive  relation  at  the  composite 
level)  using  the  material  information  at  the  subcell  level. 
To  this  end,  stress  and  strain  compatibility  equations  are 
applied.  First  of  all,  the  stress  compatibility  is 
considered  at  the  subcell  interfaces  as  given  below. 
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Strain  compatibility  is  assumed  to  be 


£f,  +Sn  =£',  +sf 
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Substitution  of  Equation  (15)  into  the  equations  in  (19) 
yields  the  stress  compatibility  equations  in  terms  of 
subcell  strains.  Therefore,  the  equations  in  (19)  become 
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These  equations  along  with  Equations  (17)  and  (20)  can 
relate  the  subcell  strains  to  unit  cell  strains-  That  is. 


M{s}  =  {s} 


(22) 
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in  which  {s}  is  a  vector  consisting  of  subcell  strains  s“ 

while  {s}  is  a  vector  consisting  of  the  unit  cell  strains  sy 
and  zeros.  Solving  Equation  (22)  results  in  the  expression 
for  the  unit  cell  strains  in  terms  of  subcell  strains.  That 
is. 


{s}  =  [^r>{s}  (23) 

Equation  (23)  gives  the  explicit  relationship  between  the 
subcell  and  the  unit  cell  strains  (i.e.  strains  at  the 
constituent  level,  the  particle  and  the  matrix,  and  strains 
at  the  composite  level.)  Equation  (23)  is  very  useful  for 
damage  and  failure  evaluation  at  the  constituent  level  as 
will  be  discussed  in  the  later  chapter.  Finally,  in  order 
to  determine  the  constitutive  relation  between  the  unit  cell 
stresses  and  strains.  Equation  (15)  is  substituted  into 
Equation  (16).  This  operation  yields  the  expression  for  the 
unit  cell  stresses  in  terms  of  the  subcell  strains.  The 
subcell  strains  are  replaced  by  unit  cell  strains  using 
Equation  (23).  AS  a  result,  the  unit  cell  stresses  become 
directly  related  to  the  unit  cell  strains.  The  relating 
constitutive  expression  is  given  in  terms  of  subcell 
material  properties  Ey,„  and  subcell  volume  fractions  7“  . 
In  other  words,  the  composite  level  constitutive  equation  is 
expressed  in  terms  of  the  particle  and  matrix  properties  as 
given  below. 


a  ij  =  EijkiZu 


(24) 
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where  Eyki  is  the  constitutive  relation  of  the  composite. 


C .  DAMAGE  MECHANICS 

Application  of  the  damage  mechanics  to  the  constituent 
materials,  the  matrix  and  the  particle  inclusions,  is 
necessary  in  order  to  predict  the  damage  initiation  and 
evolution  until  the  point  of  failure.  By  this  application, 
the  damage  modes  can  be  prescribed  fundamentally  and 
discriminately  as  matrix  cracking,  particle  cracking  or 
particle/matrix  interface  debonding .  The  present  study 
considers  matrix  cracking  with  particle  inclusions  that  are 
almost  rigid. 

For  the  development  of  the  damage  mechanics  theory,  the 
matrix  material  is  considered.  In  the  development,  for  the 
sake  of  simplicity,  vector  and  matrix  notations  will  be 
adopted  instead  of  tensor  notations.  Therefore,  stresses 
and  strains  are  expressed  in  vector  form,  and  the  strain 
increment  in  the  material  may  be  written  as, 

{dz}  =  {de}^+{ds},  (25) 

where  {ds}  is  the  total  strain  increment,  and  the  subscripts 
V  and  d  denote  the  natural  and  damaged  states,  respectively. 
In  short,  the  total  strain  increment  consists  of  the  strain 
increment  without  damage  and  the  strain  increment  as  a 
result  of  damage.  The  constitutive  equation  for  the 
undamaged  material  is  as  follows : 

{dc!)  =  [E],{dz).  (26) 
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while  the  incremental  strain  due  to  the  damage  is  assumed  to 
be  expressed  as 


{ds}. 


.. 

dA-  ‘ 


d{<T} 


(27) 


in  which  is  the  damage  potential  function.  The  damage 
potential  function  varies  depending  on  the  amount  of  damage 
state.  Thus,  it  can  be  assumed  that 


(28) 

If  the  damage  potential  function  is  also  used  for  the  damage 
criteria,  the  damage  state  must  satisfy  F^=0.  As  a  result, 
the  Equation  (28)  can  be  written  as 

^{<7}  ^{s} 

Using  Equations  (25),  (26)  and  (29),  one  can  solve  for  dX  as 
follows : 


dA  = 


{QVimds} 
{RY{Q}  +  {QV{EUQ} 


(30) 


where 


and 


{0={^} 

d{o-} 


(31) 


36 


(32) 


{R}  =  { 


a{s} 


} 


Substitution  of  Equation  (30)  into  Equation  (27)  and 
combining  the  resultant  equation  with  Equations  (25)  and 
(26)  gives 


{do}  =  [E]Ade} 


(33) 


in  which  the  constitutive  matrix  for  the  damaged  material  is 


[Eh  =  [EU[r\- 


{Q}{QV[RI  s 
{RV{Q}  +  {QV[euQ} 


(34) 


The  damage  potential  function  is  dependent  on  the 
constituent  material  and  can  be  rewritten  as 


F,  =  f({(y})-gi(zAe) 


(35) 


in  which  f  and  g  are  assumed  to  be 


(36) 


and 


g  =  (Trf +aitanh(a2(s^)J  (37) 

Here,  a ^  is  the  damage  initiation  strength  of  the  material, 
(sAe  effective  damage  strain,  and  and  atj 
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VI.  MACROMECHANICS  RESULTS  AND  DISCUSSION 


Damage  initiation,  growth,  and  saturation  were  observed 
from  the  analysis  as  the  load  increased.  Again,  Figure  6 
shows  the  stress-strain  curve  of  the  material  used  in  the 
study  and  the  finite  element  mesh  is  shown  in  Figure  3.  In 
order  to  simulate  the  experimental  condition,  a  continuous, 
uniform  strain  was  applied  to  the  specimen  in  the  modeling. 
In  examining  the  circular  hole,  it  was  noted  that  the 
saturated  damage  zone  at  a  lower  value  of  the  applied  load 
was  an  elliptical  shape  with  the  major  axis  parallel  to  the 
load  axis.  Initially,  the  major  diameter  was  much  greater 
than  the  minor  diameter.  As  the  damage  evolved,  however, 
the  minor  diameter  grew  faster  than  the  major  diameter  and 
the  damage  zone  spread  both  in  the  load  and  crack 
directions.  Figures  10,  11  and  12  shows  this  effect 
graphically.  Figure  10  denotes  the  two  regions  in  which  the 
contour  lines  of  the  crack  initiation  and  propagation  are 
plotted.  As  can  be  seen  from  the  figure.  Region  I  denotes 
the  area  up  to  the  eighth  element  in  the  x  and  y-directions 
which  is  0.162  inches  measured  from  the  center  of  the  hole. 
Region  II  denotes  the  area  up  to  the  twelfth  element  in  the 
X  and  y-directions  which  is  0.308  inches  from  the  center  of 
the  hole.  Figure  11  is  the  magnified  area  of  the  Region  I 
with  section  A  denoting  the  area  in  which  the  crack  has 
initiated.  From  the  shape  of  section  A,  it  seems  as  though 
the  major  axis  of  the  crack  is  in  the  y-direction  as  opposed 
to  in  the  x-direction  which  is  contrary  to  the  statement 
earlier.  However,  in  looking  at  the  initial  crack  area  in 
relation  to  the  area  of  the  entire  model,  it  is  a  very 
minute  section.  Therefore,  Section  A  in  Figure  11  is  simply 
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Region  I 


Region  II 


Figure  10.  Regions  Magnified  in  Figures  11 
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above  192  psi 
185  -  192  psi 
180  -  184  psi 
175  -  179  psi 
165  -  174  psi 
below  165  psi 


A  = 
B  = 
C  = 
D  = 
E  = 
F  = 


Figure  11.  Magnified  Region  I  with  Contour  Lines  of  Various 

Stress  Levels 


Figure  12.  Magnified  Region  II  with  Contour  Lines  of  Various 

Stress  Levels 
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showing  the  location  of  the  crack  formation.  It  is  not 
intended  to  show  the  shape  of  the  crack  propagation.  To 
view  the  actual  shape  and  the  effect  of  the  crack 
propagating  in  both  vertical  and  horizontal  direction, 
Figure  12  is  provided.  Figure  12  is  the  magnified  area  of 
Region  II,  and  it  shows  section  A  whose  area  represent  the 
propagation  of  the  crack,  and  as  can  be  seen,  it  is  in  both 
the  horizontal  and  the  vertical  directions.  From  the 
information  obtained  through  the  numerical  analysis 
concerning  the  final  crack  size  at  the  load  of  .2  in/in 
strain,  the  experimental  crack  size  was  compared.  As  seen 
in  Figure  13,  the  specimen  compared  well  between  the 
experimental  and  the  numerical  results.  The  experiment 
showed  cracks  of  .795  cm  (.313  in.)  and  .953  cm  (.375  in.) 
at  the  two  sides  of  the  hole  while  the  numerical  simulation 
resulted  in  crack  of  .820  cm  (.323  in.).  The  predicted 
crack  was  between  the  two  experimental  crack  sizes.  In  the 
simulation,  a  quarter  of  the  specimen  was  modeled  as  shown 
in  Figure  3  due  to  syiranetry.  However,  as  shown  by  the 
slight  difference  in  the  crack  sizes  between  the  two  sides, 
the  actual  specimen  was  not  quite  symmetric.  There  are  many 
possible  reasons  for  this  asymmetry  such  as  the  location  of 
the  hole  within  the  nonuniform  material,  the  shape  of  the 
hole,  and  the  nonunifrom  surface  roughness  of  the  hole  edge. 
A  nonuniform  particle  distribution  could  be  another  factor 
causing  such  an  asymmetry. 

Figure  14  illustrates  the  simulated  crack  propagation 
as  a  function  of  the  applied  load  for  a  specimen  with  an 
1.27  cm  (.5  in.)  diameter  hole.  The  simulation  indicated 
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Predicted  Crack  Size  =  .340  In. 


0.375  in.  0.313  in. 


Figure  13.  0.25  in.  Diameter  Hole.  Experimental  with 

Predicted  Results 


44 


that  the  resultant  crack  size  increased  with  the  increased 
applied  load.  Until  the  crack  size  became  about  10  percent 
of  the  hole  radius,  which  will  be  referred  to  as  the 
critical  crack  size  in  the  following  discussions,  the 
relationship  between  the  applied  load  and  the  crack  size  was 
bilinear.  At  the  critical  cracks  size,  however,  the 
increased  loading  propagated  and  accumulated  the  damage  zone 
in  the  direction  other  than  the  crack  orientation,  mostly  in 
the  direction  normal  to  the  crack.  Once  there  was  a 
substantial  increase  in  the  load,  the  crack  resumed 
propagation  beyond  the  critical  size  in  the  direction  of  the 
crack  orientation.  Throughout  the  region,  the  crack 
propagation  was  linearly  proportional  to  the  applied  load. 
The  slope  of  the  load  verses  the  cracks  size  was  less  after 
the  critical  crack  than  the  slope  prior  to  the  critical 
crack.  The  temporary  halt  of  the  crack  propagation,  shown 
as  the  vertical  line  segment  in  the  graph,  was  an  indication 
of  the  crack  tip  blunting  process.  As  mentioned  in  the 
earlier  section,  this  phenomena  was  also  observed  in  the 
experimental  study.  (Refer  to  the  second  figure  in  Figure 
5 . )  Recall  how  the  damage  zone  concentrated  around  the 
crack  tip  causing  the  crack  to  become  blunt  as  shown  in  the 
figure  by  a  gray,  circular  patch  of  dots  at  the  end  of  the 
crack  tip. 

Figure  15  is  the  tensile  test  coupon  that  was  used  to 
generate  the  stress-strain  curve  in  Figure  6.  The  figure 
shows  that  the  strain  experienced  is  not  uniform  within  the 
specimen.  As  shown,  the  specimen  is  divided  into 
subsections  A  through  E,  positioned  arbitrarily  along  the 
length  of  the  coupon.  The  specimen  was  then  stretched  up  to 
the  strain  value  of  .2  in/in.  Figure  15  shows  couple  of 
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Figure  15.  Deformation  of  a  Tensile  Specimen 
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representative  pictures  of  how  the  specimen  looked  as  the 
specimen  under  went  loading.  The  next  figure.  Figure  16 
shows  the  strain  verses  time  at  each  of  the  subsections 
profiled  against  the  strain  of  the  whole  length.  As  can  be 
seen,  the  deformations  at  the  various  regions  differed 
significantly.  One  of  the  major  causes  of  the  nonuniform 
deformations  is  the  nonuniform  particle  distribution  in  the 
composite.  The  largest  deformation  occurs  in  the  section  CD 
where  the  curve  is  the  steepest  and  the  highest.  The  lowest 
position  of  the  curve  DE  in  the  graph  indicate  the  least 
amount  of  deformation.  This  can  be  interpreted  as  the 
region  having  most  amount  of  particle  inclusions.  The 
overall  deformation  is  indicated  by  the  curve  AE.  As 
expected,  the  overall  curve  is  nested  between  the  subsection 
curves.  Finally,  all  the  curves  in  the  graph  follows  a 
linear  path  in  the  deformation  process  except  curve  CD. 
Curve  CD  initially  follows  a  path  that  has  a  slope  similar 
to  those  of  other  curves  until  about  two  seconds  into  the 
deformation  process.  All  of  a  sudden,  the  slope  sharply 
increases  indicating  a  greater  amount  of  deformation  per 


given  time.  As  can  be  seen  in  Figure  11,  section  CD 
represents  the  small  section  in  the  middle  of  the  sample. 
When  the  load  is  applied,  this  is  the  section  where  damage 
usually  occurs.  As  a  result,  after  a  point  in  which  the 
damage  develops,  the  deformation  of  the  material  becomes 
more  rapid  and  greater  in  magnitude. 

Clearly,  due  to  the  behavior  observed  in  the  sample 
material,  a  nonuniform  property  should  be  assumed  in  the 
modeling.  Therefore,  for  the  specimen  with  a  .635  cm  (.25 
in.)  diameter  hole,  this  was  the  case.  In  the  next  phase  of 
the  study,  it  was  hypothesized  that  although  the  overall 
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stress-strain  curve  remains  the  same,  there  can  exist  many, 
different  variations  of  particle  distribution  all  having  the 
same  overall  stress-strain  curve.  In  the  modeling,  two 
different  variations  of  65  percent  particle  volume  fraction 
(PVF)  were  used.  The  subdivisions  of  the  model  consisted  of 
two  regions  divided  along  the  line  parallel  to  the 
horizontal  axis  whose  areas  are  approximately  the  same. 
Specifically,  the  first  variation  consisted  of  .55  particle 
volume  fraction  in  the  lower  half  (the  half  where  the  hole 
resides)  with  .75  particle  volume  fraction  in  the  upper 
half.  In  the  second  variation,  the  .75  particle  volume 
fraction  was  placed  in  the  lower  half .  The  overall  particle 
volume  fraction  for  both  of  these  variations  were  .65  PVF. 
The  overall  stress-strain  curve  for  each  PVF  was  estimated 
using  a  micromechanical  model  and  the  graphs  are  shown  in 
Figures  17  and  18.  The  crack  growth  as  a  function  of 
applied  load  for  a  perforated  specimen  was  compared  between 
the  two  cases.  This  is  shown  in  Figure  19.  In  the  figure, 
"low"  denotes  the  lower  half  of  the  domain  including  the 
hole  while  "high"  denotes  the  upper  half  of  the  domain,  the 
half  where  the  load  is  applied.  The  results  indicated  that 
the  crack  behavior  around  the  notch  tip  strongly  depended  on 
the  material  property  around  the  notch,  more  specifically, 
the  particle  distribution.  Moreover,  the  graph  shows  that 
for  the  case  where  the  lower  half  is  .75  PVF,  the  strain 
(load)  reaches  up  to  the  maximum  value  of  .2  in/in,  and  the 
crack  has  advanced  as  far  as  it  will  go  at  a  crack  size/hole 
radius  value  of  .75.  The  lower  graph,  the  one  with  .55  PVF 
in  the  lower  half,  shows  that  the  crack  has  propagated 
further  along  than  the  upper  graph  under  much  less  load. 
For  the  purpose  of  looking  at  the  initial  crack  propagation. 
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the  graph  was  terminated  at  around  the  point  where  the  crack 
reached  14  elements  along  the  bottom  edge  of  the  model.  If 
we  were  to  look  at  the  whole  model,  however,  the  second 
graph  extends  all  the  way  across  indicating  complete  failure 
in  the  model  before  reaching  the  maximum  load  of  .2  in/in 
strain.  In  Figure  20,  the  crack  propagation  for  the  two 
volume  fractions  .55  and  .75  are  graphed.  Instead  of  having 
mixed  volume  fractions,  these  are  the  cases  when  the  entire 
model  was  homogeneously  at  .55  PVF  and  .75  PVF.  As  can  be 
seen,  the  two  curves  terminate  before  reaching  the  maximum 
load  of  .2  in/in  strain.  If  the  curves  are  extended,  just 
as  in  the  mixed  case  where  .55  PVF  was  in  the  lower  half, 
the  model  undergoes  a  complete  failure.  Why  is  it  then  that 
when  the  volume  fraction  is  homogenous  at  .55  PVF  or  .75  PVF 
the  model  fails  completely  but  when  they  are  mixed  together, 
(specifically,  .75  PVF  in  the  lower  half  and  .55  PVF  in  the 
upper  half)  crack  only  propagates  to  a  point  significantly 
before  the  complete  failure?  This  is  because  two  competing 
factors  are  at  work.  In  the  case  of  .55  PVF,  the  material 
is  relatively  soft,  so  it  experiences  less  stress  for  the 
same  deformation,  but  the  strength  of  the  material  is  low  as 
well.  Therefore,  the  load  overcomes  the  strength  and 
failure  occurs.  Similarly,  in  the  case  of  .75  PVF,  since 
the  material  is  stiffer,  for  the  same  amount  of  load  the 
stress  experienced  is  much  higher.  However,  the  strength  of 
the  material  is  higher  as  well  but  not  enough  to  overcome 
the  load,  and  so,  the  failure  occurs.  When  the  volume 
fractions  are  mixed  in  the  combination  of  .75  PVF  in  the 
lower  half  and  .55  PVF  in  the  upper  half,  the  stress  it 
experiences  is  greater  than  that  of  the  homogenous  case  with 
volume  fraction  of  .55  but  lower  than  that  of  .75  PVF  which 
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is  as  expected.  But  this  time,  the  stress  has  decreased 
enough  so  that  the  strengthening  effect  of  the  .75  PVF 
overcomes  the  load.  This  results  in  the  crack  propagation 
being  arrested  before  reaching  complete  failure. 

The  next  study  investigated  the  interaction  of  a  notch 
tip  and  a  pre-existing  local  defect  such  as  a  local  cluster 
of  particles  or  an  inclusion  of  a  void.  Therefore,  a  local 
zone  ahead  of  the  notch  tip  was  assumed  to  have  a  much 
higher  stiffness  or  a  much  lower  stiffness  that  the 
composite  material  depending  on  the  type  of  inclusion.  They 
are  represented  by  the  shaped  regions  in  the  figure. 
Specifically,  this  was  modeled  with  the  stiffness  being 

greater  or  lesser  by  a  factor  of  10® .  These  effects  were 
examined  on  the  initiation  and  growth  of  a  crack  from  the 
notch  tip.  Figure  21  shows  the  locations  of  the  inclusion. 
Specifically,  they  were  placed  at  .183  cm  (.072  in.),  .411 
cm  (.162  in.),  and  .953  cm  (.375  in.)  from  the  edge  of  the 
hole.  Figure  22  is  a  plot  of  the  crack  growth  with  a  void 
at  two  different  locations .  When  the  void  was  located  at 
the  distance  of  2.75  times  the  hole  radius  (i.e.  void  at 
.345  in.),  the  void  did  not  affect  the  crack  initiation  and 
growth  at  the  early  stage  because  of  the  distance  being  too 
far  from  the  place  of  crack  initiation.  However,  as  the 
crack  propagated  toward  the  void  and  became  closer,  the  void 
affected  the  crack  growth  more  significantly.  On  the  other 
hand,  when  the  void  was  placed  at  a  distance  much  closer  to 
the  hole  edge,  at  a  distance  of  .183  cm  (.072  in.),  the 
crack  initiation  developed  at  a  much  smaller  load  (i.e. 
applied  strain),  and  the  crack  propagated  at  a  load  lower 
than  that  of  the  case  where  there  was  no  void.  The  study 
indicated  that  for  the  present  particulate  composite 


56 


Figure  21.  Locations  of  Void  or  Particle  Cluster 

Inclusions 
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Applied  Strain 


Figure  22.  Effects  of  a  Particle  Cluster  Near  Circular  Notch 


Applied  Strain 


Figure  23.  Effects  of  a  Local  Void  Near  a  Circular  Notch 


material,  an  inclusion  of  a  void  at  a  distance  about  the 
size  of  the  hole  diameter  from  the  edge  made  almost  no 
difference  in  the  initial  crack.  Figure  23  compares  the 
crack  growth  between  the  cases  with  a  void  and  a  particle 
cluster,  respectively,  at  the  same  distance  from  the  hole 
edge.  As  expected,  a  particle  cluster  prevented  the  crack 
initiation  and  its  growth  compared  to  the  no  inclusion  case. 
On  the  other  hand,  the  presence  of  a  void  assisted  in  the 
crack  initiation  and  growth. 

The  final  study  compared  the  damage  and  crack  between 
two  cases  for  a  model  with  a  hole  diameter  of  1.27  cm  (.5 
in.).  The  first  case  involved  applying  the  load  with  a 
uniform  traction  which  will  be  called  traction  control,  and 
the  second  case  involved  applying  the  load  with  a  uniform 
displacement  load  which  will  be  called  displacement  control. 
The  two  cases,  although  different  in  form,  had  the  same 
amount  of  total  applied  load.  Figure  24  shows  the 
comparison  of  the  boundary  edge  traction  between  the 
traction  control  and  the  displacement  control  cases .  The 
displacement  control  begins  with  a  stress  value  lower  than 
that  of  the  traction  control  case  but  then  increases  as  the 
distance  along  the  x-axis  increases.  The  traction  control 
in  the  mean  while  stays  constant.  In  Figure  25,  the 
comparison  of  boundary  edge  displacement  is  shown.  Here, 
the  displacement  for  the  traction  control  begins  above  the 
displacement  control  then  decreases  to  a  final  value  under 
the  constant  displacement  control  curve.  These  two  figures 
collectively  show  that  the  traction  control  resulted  in  a 
larger  load  at  the  hole  location,  and  also,  the  boundary 
edge  displacement  was  greater  as  well.  As  a  result,  the 
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traction  control  yielded  a  greater  deimage  zone  and  a  longer 
crack  size  in  the  vicinity  of  the  hole. 
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Displacement  (in) 


Figure  25.  Comparison  of  Boundary  Edge  Displacement  for  a 
Performated  Plate  Between  the  Traction  Control  and  the 

Displacement  Control 


VII.  MICRO /MACROMECHANICS  RESULTS  AND  DISCUSSION 


The  analysis  used  in  this  approach  simulated  the  matrix 
damage  initiation  and  damage  propagation  using  the 
micromechanical  model  and  the  damage  mechanics  as  described 
in  Chapter  IV.  When  the  damage  was  saturated  in  a  local 
zone,  the  location  could  not  sustain  any  additional  loading. 
Therefore,  it  was  assumed  that  a  crack  occurred  within  the 
saturated  damage  zone.  In  other  words,  the  crack  length  was 
deemed  to  be  equal  to  the  size  of  the  saturated  damage  zone. 
Again,  the  stress-strain  curve  for  the  material  used  in  the 
study  is  provided  in  Figure  26.  This  curve  represents  the 
macro-level,  smeared  material  property.  The  micromechanical 
model  requires  the  material  properties  of  the  constituent 
materials  and  their  volume  fractions.  However,  initially, 
this  information  is  not  known.  Therefore,  it  is  necessary 
to  back  calculate  the  micro-level  material  properties  from 
the  macro-level  material  properties  using  the  proposed 
micromechanical  model  and  the  damage  mechanics .  The 
computed  micro-level  material  properties  were  as  follows: 

•  particle  elastic  modulus  (10  x  10^  psi) 

•  matrix  elastic  modulus  (335  psi) 

•  volume  fraction  (.65) 

In  addition,  the  damage  initiation  strength  of  the  matrix 
material  calculated  using  equation  (37)  in  Chapter  IV  was  70 
psi,  with  damage  material  constants  being  a,  =92  psi  and 
^2=12  psi-  The  reproduced  stress-strain  curve  using  the 
micromechanical  model  and  the  damage  mechanics  is 
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superimposed  on  the  experimental  stress-strain  curve  in 
Figure  26.  As  can  be  seen,  the  experimental  curve  and  the 
predicted  curve  corresponded  well  to  one  another.  By  this 
comparison,  the  micro-level  material  values  were  deemed  to 
be  reasonably  accurate. 

The  composite  material  was  now  tested  with  a  hole  in 
the  center.  Tensile  loads  were  applied  to  the  specimens 
with  constant  crosshead  speed  of  .2  in/min.  Cracks 
initiated  from  stress  concentration  at  the  edges  of  the 
holes.  Figures  27  through  29  show  the  cracks  in  the 
specimens  at  strain  of  .2  in/in.  Two  of  the  specimens  had 
.25  in.  diameter  holes  while  the  third  specimen  had  a  .50 
in.  diameter  hole.  The  first  two  specimens  with  the  same 
diameter  hole  show  two  very  different  crack  lengths  under 
the  same  load  as  shown  in  Figures  27  and  28.  In  the  third 
specimen  (the  one  with  .50  in.  diameter  hole),  the  crack 
lengths  on  the  two  sides  of  the  hole  came  out  to  be  quite 
asymmetric,  as  shown  in  Figure  29.  These  test  results  can 
be  explained  from  the  nonhomogeneity  of  particle 
distribution  in  the  matrix  material.  As  observed  from  the 
tensile  test  of  a  dogbone  shaped  specimen  (without  a  hole), 
the  local  strain  was  very  different  from  location  to 
location  within  the  specimen  under  the  same,  constant  load. 
(Refer  to  Figure  16)  The  experimental  stress-strain  curve 
in  Figure  26  is  an  average  strain  between  the  gauge  length 
of  2.7  in.  Therefore,  when  a  hole  is  drilled  into  a 
specimen,  the  hole  may  be  drilled  into  a  particle-dense 
zone  or  into  a  particle-sparse  zone  without  affecting  the 
overall  stress-strain  curve.  Because  the  specimen  was  under 
strain  control,  the  local  strain  and  stress  around  the  hole 
varied  depending  on  what  type  of  topography  (particle-dense 
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Figure  27.  Comparison  of  Experimental  and  Theoretical  Crack 
Sizes  for  a  0.25  in.  Diameter  Hole  Sample 
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Figure  28.  Comparison  of  the  Experimental  and  Theoretical 
Crack  Sizes  for  Another  0.25  in  Diameter  Hole  Sample 
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Figure  29.  Comparison  of  Experimental  and  Theoretical  Crack 
Sizes  for  a  0.5  in.  Diameter  Hole  Sample 


or  particle-sparse)  the  hole  was  placed.  Therefore,  it  was 
quite  possible  to  have  two  different  crack  sizes  for  the 
same  material,  under  similar  conditions  which  was  the  case 
in  our  experiment.  Similarly,  the  asymmetric  cracks 
observed  in  the  third  specimen  can  be  explained  by  the  same 
reasoning  as  well  as  other  factors  such  as  asymmetric  hole 
shape,  location,  and  misaligned  loading. 

Computer  simulation  was  conducted  next  to  examine  if 
similar  results  can  be  obtained  and  thereby,  test  the 
accuracy  of  the  micro/macromechanics  method  in  predicting 
the  crack  initiation  and  propagation.  Initially,  it  was 
assumed  that  the  specimen  had  a  uniform  material  property 
(i.e.  uniform  particle  distribution).  The  specimen  with  a 
.25  in.  diameter  hole  was  considered  first.  A  quarter  of 
the  specimen  was  modeled  as  illustrated  in  Figure  3 .  The 
predicted  crack  size  obtained  from  the  computer  simulation 
is  shown  in  Figure  27  along  with  the  actual  cracks  obtained 
from  the  experiment.  The  predicted  crack  length  came  out  to 
be  .340  in.  while  the  experimental  crack  lengths  were  .375 
in.  and  .313  in.  The  predicted  crack  length  corresponded 
well  to  the  experimental  values  for  this  sample. 

In  the  simulation  for  the  second  specimen,  a  mixture  of 
particle  volume  fractions  was  used.  Figure  30  shows  the 
division  of  the  model.  Two  cases  were  considered.  The 
first  case  considered  the  lower  half  of  the  specimen  at  a 
PVF  of  .55  and  the  upper  half  at  a  PVF  of  .75.  In  the 
second  case,  the  particle  volume  fractions  were  switched  so 
that  the  lower  half  had  a  PVF  of  .75  while  the  upper  half 
had  a  .55  PVF.  However,  in  both  cases,  the  average  PVF  of 
.65  was  maintained.  The  particle  volume  fraction  (PVF) 
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0.69 
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Figure  30.  Shows  the  Division  of  the  Model  Between  the  Upper 

Half  and  the  Lower  Half 
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represent  the  degree  of  particle  inclusion  in  the  matrix 
material  such  that  greater  the  PVF,  greater  the  amount  of 
particle  inclusion.  So,  in  the  first  case  where  the  hole 
lied  in  the  particle-sparse  zone,  the  predicted  crack  size 
came  out  to  be  .063  in.  which  exactly  matched  one  of  the 
experimental  crack  lengths.  The  other  experimental  crack 
length  was  .036  in.  On  the  other  hand,  the  second  case 
yielded  a  predicted  crack  length  of  .460  in.  which  was 
significantly  higher  than  either  of  the  experimental  values. 
This  indicates  that  the  cracks  produced  in  the  experiment 
was  due  to  the  hole  being  drilled  into  a  particle-sparse 
zone. 

Some  other  related  studies  were  undertaken  next  to 
investigate  the  crack  initiation  and  propagation  under 
different  conditions.  This  time,  different  particle  volume 
fractions  were  considered.  The  particle  distribution  was 
assumed  to  be  uniform  across  the  specimen  with  a  .25  in. 
diameter  hole  at  various  PVF's.  These  models  were  then 
placed  under  the  same  loading  conditions  as  before,  and  the 
results  were  obtained.  The  crack  size  verses  applied  strain 
is  plotted  in  Figure  31.  As  expected,  the  graphs  shows  that 
the  higher  the  PVF,  the  larger  the  crack  size  under  the 
sample  applied  strain.  One  thing  to  be  noted  here  is  that 
when  the  PVF  was  .7  or  higher,  the  crack  size  increased 
continuously  without  a  halt  as  the  specimen  was  stretched. 
On  the  other  hand,  when  the  PVF  was  .65  or  lower,  there  was 
a  temporary  halt  in  the  crack  propagation  as  the  specimen 
was  stretched.  Specifically,  the  crack  did  not  propagate  at 
certain  lengths  until  there  was  a  significant  jump  in  the 
applied  load.  This  phenomena  was  interpreted  to  indicate  a 
crack  blunting  process.  Also,  when  the  PVF  was  .55,  a  large 
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Applied  Strain 


Figure  31.  Crack  Propagation  of  Various  Volume  Fractions 
Under  Displacement  Loading 
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increase  in  the  applied  strain  was  needed  to  initiate  a 
crack. 

The  next  study  conducted  was  the  response  due  to  stress 
control.  Uniform  traction  was  applied  at  the  top  edge  of 
the  specimen  and  was  increased  continuously.  The  crack 
growth  with  different  PVF's  is  given  in  Figure  32.  The 
results  indicated  that  the  crack  size  increased  continuously 
at  almost  the  same  rate  as  the  applied  load  regardless  of 
PVF.  Under  this  type  of  loading,  the  stress  state  was 
similar  that  there  was  not  much  difference  in  the  crack  size 
among  the  various  PVF's.  A  higher  PVF  yielded  a  lower 
stress  in  the  matrix,  because  the  particles  supported  a 
larger  portion  of  the  load.  Hence,  the  higher  PVF,  the 
smaller  the  crack  size.  But  again,  the  difference  was 
hardly  noticeable  and  was  much  less  than  that  of  the 
displacement  control  models.  This  outcome  was  as  predicted. 
Since  the  load  was  in  the  form  of  a  uniform  stress,  the 
stress  experienced  throughout  the  material  was  uniform 
although  the  deformation  changed  from  point  to  point. 
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Figure  32.  Crack  Propagation  of  Various  Volume 
Fractions  Under  Traction  Loading 


VIII.  CONCLUSIONS  AND  RECOMMENDATIONS 


Modeling  of  crack  initiation  and  propagation  in  a 
particulate  composite  (solid  rocket  propellant)  material  was 
conducted  using  macromechanical  and  micro/macromechanical 
theories.  The  macromechanical  studies  showed  that  the  crack 
behavior  obtained  from  the  experiment  could  be  reproduced 
using  finite  element  methods  and  a  smeared,  homogenous 
material  property  representing  the  particulate  composite. 
Also,  the  interpretation  of  damage  growth  and  saturation 
observed  in  the  model  as  being  the  initiation  and 
propagation  of  a  crack  was  accurate.  In  the 
micro/macromechanical  model  and  the  damage  mechanics,  it  was 
shown  that  the  micro/macromechanical  approach  in  predicting 
matrix  cracking  in  particulate  composites  was  appropriate 
and  accurate.  Furthermore,  this  approach  was  able  to  easily 
model  the  nonuniformity  in  particle  distribution  in  the 
composite  material.  The  predication  of  the  crack  size  in  a 
nonuniformly  distributed  composite  was  in  good 
correspondence  with  the  experimental  results.  Further 
analysis  indicated  that  the  displacement  control  of  tensile 
specimens  resulted  in  a  smaller  crack  size  for  a  smaller  PVF 
subjected  to  the  same  global  strain.  On  the  other  hand,  the 
load  controlled  specimens  yielded  quite  the  opposite 
results.  However,  the  displacement  controlled  specimens 
produced  a  much  larger  difference  among  the  different  PVF's 
than  the  load  controlled  specimens. 

Recommendations  for  further  study  are  as  follows; 
Incorporate  time  depnedent  effects  such  as  strain  rate 
effect  and  viscous  effects,  include  more  complex  load 
conditions  such  as  biaxial  and  shear  loading,  and  conduct 
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analysis  on  an  actual,  operational  propellant  material.  The 
composite  used  in  this  study  was  an  inert  material  with  a 
similar  composite,  material  structure. 
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APPENDIX 


SAMPLE  INPUT  FOR  TRACTION  LOAD  (MICROMECHANICS) 

806  1  17  2  0  2  10  0  0 
0.0  19.  1.  1.  1.0e-4  1.0e-4  1.0 
0  0  0  0  0  0 

0  0  0  0 

0  750  0  0 
0  0 
0  0 
1  6 

1.  10000000.  335.  0.2  0.5  2000.  70.5000.7000.0.65 
1.  1.  90.0  12.0  0.  0.  0.  0.  0.  0. 

0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 


1 

5 

0 

0.125000E+00  .0.426851  E-17 

O.OOOOOOE+00 

2 

5 

0 

0.134329E+00  O.OOOOOOE+00 

O.OOOOOOE+00 

3 

5 

0 

0.145435E+00  O.OOOOOOE+00 

O.OOOOOOE+OO 

4 

5 

0 

0.159176E+00  O.OOOOOOE+00 

O.OOOOOOE+OO 

• 

* 

* 

*  <t 

* 

* 

* 

* 

*  * 

it 

803 

6 

0 

-0.905402E-17 

0.104009E+01 

O.OOOOOOE+OO 

804 

6 

0 

-0.627097E.17 

0.114304E+01 

O.OOOOOOE+OO 

805 

6 

0 

-0.325488E-17 

0.125460E+01 

O.OOOOOOE+OO 

806 

6 

0 

-0.213426E-15 

0.137500E+01 

O.OOOOOOE+OO 

1 

1 

806 

780  779  805 

2 

1 

805 

779  778  804 

3 

1 

804 

778  777  803 

4 

1 

803 

777  776  802 

5 

1 

802 

776  775  801 

* 

* 

*  * 

*  * 

*  * 

*  * 

* 

746 

* 

1 

*  « 

760 

*  * 

786  787  761 

747 

1 

761 

787  788  762 

748 

1 

762 

788  789  763 

749 

1 

763 

789  790  764 

750 

1 

764 

790  791  765 

390  2 
0.  0. 

20.  9.375 
416  2 
0.0 

20.  18.75 
442  2 
0.  0. 

20.  18.75 
468  2 
0.  0. 


20.  18.75 
780  2 
0.  0 

20.  18.75 
806  2 
0.  0. 

20.  9.375 


79 


SAMPLE  INPUT  FOR  DISPLACEMENT  LOAD  WITH  MULITPLE  MATERIAL 


806  2  17  2  0  1  10  0  0 
0.0  10.  1.  1.  1.0e-3  1.0e-3  1.0 
0  0  0  0  0  0 

0  0  0  0 

0  750  0  0 
0  0 
0  0 

1  4 

1700.  7.  0.5  1.0  0.  0.  0.  1003.54  1026.59  928.38 
803.25  683.72  515.17  311.65  108.37  6.  0.  0.  0.  0. 
0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

2  4 

1700.  7.  0.5  1.0  0.  0.  0.  1893.04  1941.36  1765. 
1321.98  802.28  324.42  51.63  57.32  6.  0.  0.  0.  0. 


0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

1  5  0  0.250000E+00  -0.426851  E-1 7 

2  5  0  0.262500E+00  O.OOOOOOE+00 

3  5  0  0.275000E+00  O.OOOOOOE+00 

4  5  0  0.287500E+00  O.OOOOOOE+00 

5  5  0  0.300000E+00  O.OOOOOOE+00 

6  5  0  0.312500E+00  O.OOOOOOE+00 

7  5  0  0.325000E+00  O.OOOOOOE+00 


O.OOOOOOE+00 

O.OOOOOOE+00 

O.OOOOOOE+00 

O.OOOOOOE+00 

O.OOOOOOE+00 

O.OOOOOOE+00 

O.OOOOOOE+00 


803  6  0  0.157596E-16  0.104009E+01 

804  6  0  0.109154E-16  0.114304E+01 

805  6  0  0.566550E-17  0.125460E+01 

806  6  0  -0.213426E-15  0.137500E+01 

1  1  806  780  779  805 

2  1  805  779  778  804 

3  1  804  778  777  803 

4  1  803  777  776  802 

5  1  802  776  775  801 

6  1  801  775  774  800 


O.OOOOOOE+00 

O.OOOOOOE+00 

O.OOOOOOE+00 

O.OOOOOOE+00 


748  1 

749  1 

750  1 
390  2 

0.  0. 

10.  2.948 
416  2 
0.  0. 

10.  5.896 
0.  0. 


* 

* 

* 

He 

* 

* 

762 

788 

789 

763 

763 

789 

790 

764 

764 

790 

791 

765 

780  2 
0.  0 

10.  5.896 
806  2 
0.  0. 

10.  2.948 
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SAMPLE  INPUT  FOR  DISPLACEMENT  LOAD  (MICROMECHANICS) 


806  1  0  0  0  0  1  0  17 
0.0  100.  1.  1.  1.0e-4  1.0e-4  1.0 
0  0  0  0  0  0 
0  0  0  0 
0  750  0  0 
0  0 
0  0 
1  6 

1.  10000000.  335.  0.2  0.5  2000.  70.5000.7000.0.65 


1. 

1.  90.0  12.0  0.  0.  0.  0. 

p 

p 

0. 

0.  0. 

0.  0, 

.  0.  0.  0.  0.  0. 

0. 

0.  0. 

0.  0, 

,  0.  0.  0.  0.  0. 

1 

5 

0 

0.250000E+00 

-0.426851  E-1 7 

O.OOOOOOE+00 

2 

5 

0 

0.262500E+00 

O.OOOOOOE+00 

O.OOOOOOE+00 

3 

5 

0 

0.275000E+00 

O.OOOOOOE+00 

O.OOOOOOE+00 

4 

5 

0 

0.287500E+00 

O.OOOOOOE+00 

O.OOOOOOE+OO 

5 

5 

0 

0.300000E+00 

O.OOOOOOE+00 

O.OOOOOOE+00 

6 

5 

0 

0.312500E+00 

O.OOOOOOE+00 

O.OOOOOOE+OO 

7 

* 

5 

* 

0 

* 

0.325000E+00 
*  * 

O.OOOOOOE+00 

* 

O.OOOOOOE+OO 

* 

* 

* 

*  * 

* 

* 

* 

* 

*  * 

* 

801 

6 

0 

0.242738E-16 

0.8591 52E+00 

O.OOOOOOE+OO 

802 

6 

0 

0.202087E-16 

0.945540E+00 

O.OOOOOOE+OO 

803 

6 

0 

0.157596E-16 

0.104009E+01 

O.OOOOOOE+OO 

804 

6 

0 

0.109154E-16 

0.114304E+01 

O.OOOOOOE+OO 

805 

6 

0 

0.566550E-17 

0.125460E+01 

O.OOOOOOE+OO 

806 

7 

0 

-0.213426E-15 

0.137500E+01 

O.OOOOOOE+OO 

1 

1 

806 

780  779  805 

2 

1 

805 

779  778  804 

3 

1 

804 

778  777  803 

1 

1 

803 

777  776  802 

5 

1 

802 

776  775  801 

6 

♦ 

1 

* 

801 
*  tr 

775  774  800 
*  * 

* 

* 

* 

* 

*  * 

*  * 

*  * 

*  * 

747 

1 

761 

787  788  762 

748 

1 

762 

788  789  763 

749 

1 

763 

789  790  764 

750 

1 

764 

790  791  765 

390  2  .25 
416  2.25 
442  2  .25 


754  2  .25 
780  2  .25 
806  2  .25 
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SAMPLE  INPUT  FOR  DISPLACEMENT  LOAD  (MACROMECHANICS  CASE) 


806  2  0  0  0  2  1  0  17 
0.0  100.  1.  1.  1.0e-3  1.0e-3  1.0 
0  0  0  0  0  0 
0  0  0  0 
0  750  0  0 
0  0 
0  0 

1  4 

1700.  13.  0.5  1.0  0.  0.  0.  1294.  110.  1142.9 
130.  930.2  150.  689.7  170.  512.8  180.  220.2  192.  6. 

0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

2  1 

1000.  12940000.  0.5  1.0  0.  0.  0.  0.  0.  0. 

0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

0.  0.  0.  0.  0.  0.  0.  0.  0.  0. 

1  5  0  0.125000E+00  -0.426851  E-1 7  O.OOOOOOE+00 

2  5  0  0.134329E+00  O.OOOOOOE+00  O.OOOOOOE+00 

3  5  0  0.145435E+00  O.OOOOOOE+00  O.OOOOOOE+00 

4  5  0  0.159176E+00  O.OOOOOOE+00  O.OOOOOOE+00 

5  5  0  0.176205E+00  O.OOOOOOE+OO  O.OOOOOOE+00 

6  5  0  0.197077E+00  O.OOOOOOE+OO  O.OOOOOOE+OO 


* 

* 

* 

* 

« 

* 

* 

* 

* 

* 

* 

* 

* 

* 

803 

6 

0 

-0.905402E-17 

0.104009E+01 

804 

6 

0 

-0.627097E-17 

0.114304E+01 

805 

6 

0 

-0.325488E-17 

0.125460E+01 

806 

7 

0 

-0.213426E-15 

0.137500E+01 

1 

1 

806 

780 

779  805 

2 

1 

805 

779 

778  804 

3 

1 

804 

778 

777  803 

1 

1 

803 

777 

776  802 

5 

1 

802 

776 

775  801 

6 

1 

801 

775 

774  800 

7 

1 

800 

774 

773  799 

• 

* 

*  * 

*• 

* 

*  * 

* 

* 

* 

* 

*  * 

* 

* 

748 

1 

762 

788 

789  763 

749 

1 

763 

789 

790  764 

750 

1 

764 

790 

791  765 

390  2  .2751 
416  2.2751 
442  2  .2751 


O.OOOOOOE+OO 

O.OOOOOOE+OO 

O.OOOOOOE+OO 

O.OOOOOOE+OO 


754  2  .2751 
780  2  .2751 
806  2  .2751 
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SAMPLE  INPUT  FOR  TRACTION  LOAD  (MACROMECHANICS  CASE) 

806  1  17  2  0  2  10  0  0 
0.0  10.  1.  1.  1.0e-4  1.0e-4  1.0 
0  0  0  0  0  0 
0  0  0  0 
0  750  0  0 
0  0 
0  0 
1  6 

1.  10000000.  335.  0.2  0.5  2000.  70.5000.7000.0.65 


1. 

1.  90.0  12.0  0.  0.  0.  0.  0.  0. 

0. 

0.  0. 

0,  0 

.  0.  0.  0.  0.  0. 

0. 

0.  0. 

0.  0 

.  0.  0.  0.  0.  0. 

1 

5 

0 

0. 1 25000E+00  -0.426851  E-1 7 

O.OOOOOOE+OO 

2 

5 

0 

0. 1 34329E+00  O.OOOOOOE+00 

O.OOOOOOE+OO 

3 

5 

0 

0. 1 45435E+00  0.  OOOOOOE+00 

O.OOOOOOE+OO 

4 

5 

0 

0. 1 59 1 76E+00  O.OOOOOOE+OO 

O.OOOOOOE+OO 

5 

5 

0 

0.176205E+00  O.OOOOOOE+OO 

O.OOOOOOE+OO 

6 

5 

0 

0.1 97077E+00  O.OOOOOOE+OO 

O.OOOOOOE+OO 

7 

5 

0 

0.222285E+00  O.OOOOOOE+OO 

O.OOOOOOE+OO 

8 

5 

0 

0.252279E+00  O.OOOOOOE+OO 

O.OOOOOOE+OO 

9 

5 

0 

0.287474E+00  O.OOOOOOE+OO 

O.OOOOOOE+OO 

* 

*  *  * 

* 

* 

*  *  * 

* 

* 

* 

*  *  * 

798 

6 

0 

-0. 1 96987E-1 6  0.646339E+00 

O.OOOOOOE+OO 

799 

6 

0 

-0. 1 79823E-1 6  0.709830E+00 

O.OOOOOOE+OO 

800 

6 

0 

-0. 1 60672E-1 6  0.780671  E+00 

O.OOOOOOE+OO 

801 

6 

0 

-0.1 39455E-1 6  0.8591 52E+00 

O.OOOOOOE+OO 

802 

6 

0 

-0. 1 1 6 1 01  E-1 6  0.945540E+00 

O.OOOOOOE+OO 

803 

6 

0 

-0.905402E-1 7  0.1 04009E+01 

O.OOOOOOE+OO 

804 

6 

0 

-0.627097E-1 7  0.11 4304E+01 

O.OOOOOOE+OO 

805 

6 

0 

-0.325488E-1 7  0. 1 25460E+01 

O.OOOOOOE+OO 

806 

6 

0 

-0.21 3426E-1 5  0. 1 37500E+01 

O.OOOOOOE+OO 

1 

1 

806 

780  779  805 

2 

1 

805 

779  778  804 

3 

1 

804 

778  777  803 

4 

1 

803 

777  776  802 

5 

1 

802 

776  775  801 

6 

1 

801 

775  774  800 

7 

1 

800 

774  773  799 

* 

* 

*  * 

*  * 

* 

*  * 

*  * 

* 

* 

♦  * 

*  ♦ 

* 

* 

*  * 

*  * 

746 

1 

760 

786  787  761 

747 

1 

761 

787  788  762 

748 

1 

762 

788  789  763 

749 

1 

763 

789  790  764 

750 

1 

764 

790  791  765 

390  2 

0.  0, 

10.  3.9315 
416  2 
0.  0 

10.  7.863 


780  2 
0.0 

10.  7.863 
806  2 
0.  0. 

10.  3.9315 
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